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Abstract 


In  this  paper  we  use  combinatorial  techniques  to  solve  recurrence  relations  In 
two  variables  of  the  form 


T(N,k)  =  2  T(N/2,k)  +  T(N,k-l)  *  t(N) 

and  related  recurrences.  These  recurrences  arise  in  the  analysis  of  algorithms 
based  on  a  paradigm  called  "multidimensional  divide  -and -conquer"!  The  analyses 
that  we  present  are  interesting  from  a  combinatorial  view,  and  show  that  certain 
algorithms  are  very  efficient. 
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1.  Introduction 

In  this  paper  we  shall  study  the  problem  of  (exactly)  analyzing  Multidimensional 
Divide-and-Conquer  (MDC)  algorithms  described  by  Bentley  [1978].  The  execution  costs  of 
these  algorithms  are  usually  described  by  recurrence  relations  in  two  variables  of  the  form 

T(N,k)  -  a  T(N/2,k)  +  b  T(N,k-l)  ♦  f(N) 

or 

T(N,k)  «■  a  T(N/2,k)  ♦  b  T(N/2,k-l)  ♦  f(N) 

with  initial  values 

T(l,k)  «  (XI)  and  T(N,2)  -  g(N) 

where  a  and  b  are  integers,  and  f  and  g  are  functions  of  N.  These  recurrences  have  only 
been  roughly  solved  for  fixed  k.  The  purpose  of  this  paper  is  to  solve  these  recurrences 
exactly  using  combinatorial  techniques.  The  exact  analysis  gives  us  the  constant  factor  in  the 
expression  of  T(N,k)  as  a  function  of  k,  and  also  the  ability  to  compare  MDC  algorithms  with 
more  obvious  algorithms. 

In  Section  2  we  sketch  a  particular  MDC  algorithm  and  derive  the  recurrence  describing  its 
running  time.  We  solve  that  recurrence  precisely  in  Section  3,  and  in  Section  4  we  use  the 
combinatorial  solution  to  describe  the  behavior  of  the  algorithm.  Section  5  is  a  collection  of 
MDC  recurrences  and  their  solutions.  In  Section  6  we  extend  the  solution  method  to  more 
general  recurrences,  and  conclusions  are  offered  in  Section  7. 


2.  The  All-Points  ECOF  Algorithm 


In  this  section  we  shall  investigate  a  particular  MDC  algorithm  and  show  how  its  recurrence 
can  be  derived;  this  algorithm  is  due  to  Bentley  and  Shamos  and  is  described  in  Bentley 
[1978].  Our  purpose  in  this  section  is  not  to  learn  alt  the  details  of  the  algorithm,  but  rather 
to  understand  how  its  recurrence  arises.  We  say  that  a  point  X-fxj^x^)  in  an  Euclidian 
k-space  dominates  point  Y  iff  Xj  z  yj  for  all  i.  The  rank  r(X)  of  a  point  X  is  the  number  of 
points  dominated  by  X.  Given  N  points  in  k-space,  the  All-Points  ECDF  Problem  is  to  compute 
the  rank  of  each.  The  following  is  a  sketch  of  Algorithm  ECDFk  described  by  Bentley  [1978] 
for  solving  this  problem  on  a  set  S  of  N  points  in  k-space. 

l.If  the  number  of  points,  N,  in  S  is  one,  then  solve  the  problem  in  0(1)  time;  if  the 

dimension,  k,  of  the  point  set  is  two,  then  solve  the  problem  in  0(N  Ig  N)  time.  If 

either  of  these  conditions  holds,  return  to  the  caller;  otherwise  continue  to  Step  . 

-  Buff  Section 
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S  in  two  subsets  A  and  B,  each  containing  N/2  points. 

3.Recur$ively  solve  the  all-points  ECDF  problem  on  A  and  B  (each  problem  of  N/2 
points  in  k-space). 

^.Remarking  that  any  point  of  B  dominates  each  point  of  A  in  x-coordinate,  we  may 
remove  this  coordinate  and  solve  the  reduced  problem  (of  finding  for  each  point 
in  B  how  many  points  of  A  it  dominates)  on  the  N  points  projected  in 
(k-l)-space. 

Let  us  denote  by  T(N,k)  the  time  for  solving  the  all-points  ECDF  problem  on  a  set  of  N 
points  in  k-space.  Since  Step  2  can  be  performed  in  time  0(N),  we  find  the  recurrence 
relation 

T(N,k)-2  T(N/2,k)+T(N,k-l)+0(N).  (1) 

From  Step  1  we  have  the  boundary  conditions^ 

T(N,2)=0(N  log  N) 

T(l,k)**0(l) . 

To  analyze  the  above  recurrence  we  must  remove  the  "thetas".  We  will  therefore  solve  the 
related  recurrence 

T(N,k)-2  T(N/2,k)VT<N,k-l)*N  (2) 

T(N,2)-N  log  N 
TCI  AM- 

The  solution  of  the  Equation  2  is  always  within  a  constant  factor  of  solution  of  Equation  1, 
and  therefore  is  precise  enough  for  our  purposes. 


3.  Combinatorial  Solution 

We  shall  first  solve  Equation  2  assuming  that  N  -  2P.  For  this  purpose,  we  define 
B(p,k)-T(2p,k)/2p  . 

Dividing  both  sides  of  Equation  2  by  N  and  substituting  this  definition  of  B  gives  the  reduced 
system 

8(p,k)“B(p- 1  ,k)*B(p,k- 1  )♦  1  (3) 

B(p,2)-p 

B(0,kM. 

We  may  easily  verify  that  the  solution  to  Equation  3  is 


^Mara  Is  N  danofe*  always  lof  ^  N 
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fp+K-2\  (p«-k-3\ 
B(p,k)  -  2  \  k-2  /  ♦  V  k-1  I- 


1. 


-  3  - 

(4) 


The  solution  of  Equation  2  is  now  immediate;  it  is 

T  /lg  N  +(k-2)\  /lg  N  +(k-3)\  1 

T(N,k)  -  N  L  2  V  k-2  M  k-1  /  -  1  J. 


(5) 


We  know  now  what  the  solution  is;  we  do  not,  however,  know  why  it  is  so.  To  achieve  a 
more  intuitive  understanding  of  this  solution,  we  shall  transform  Equation  3  by  defining 
A(p,k)=B(p,k)+l  and  find  a  combinatorial  interpretation  of  the  new  recurrence 

A(p,k)-A(p-l,k)+A(p,k-l)  (6) 

A(p,2)=p+1 

A(0,k)«=2. 

Let  us  consider  the  net  of  Figure  1.  The  number  W(ab)(p,k)  of  ways  from  (a,b)  to  (p,k)  using 
only  edges  in  the  network  satisfies  the  recurrence 


W,  b(p,k)«=Wa  b(p- 1  ,k)*Wa  b(p,k- 1 ) 

W^fp.OM  ’ 

Wa;b<0,k)=l 

which  is  the  same  as  Recurrence  6,  except  for  the  bounds.  Clearly,  to  go  from  (a,b)  to  (p,k) 
we  must  choose  (p-a)  horizontal  steps  from  a  total  of  (p-a)+(k-b)  steps  (the  other  k-b  being 
vertical).  This  immediately  yields 


(p-a+k-b\ 
k-a  I. 

We  therefore  build  the  net  corresponding  to  Equation  6  and  expand  it  on  the  line  k*l{  this  is 
depicted  in  Figure  2.  It  is  now  easy  to  interpret  A(p,k)  in  terms  of  counting  paths,  and  we 
find  that 

A(p,k)  -  2  WQ  ^tp.k)  ♦  W2>i(p,k) 

which  yields  Equation  4. 

In  the  next  section,  we  shall  study  several  other  recurrences  using  this  method,  which  we 
can  summarize  as  follows. 

-Manipulate  the  initial  recurrence  until  obtaining  a  form  similar  to  Equation  6. 

-Build  the  corresponding  net,  and  (if  passible)  expand  it  to  make  the  bounds 
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Figure  1 :  A  very  simple  ret 


Figure  2:  Net  associated  with  Recurrence  6 

"more  simple". 

-Interpret  the  recurrence  in  terms  of  paths  and  solve  it. 

4.  Interpretation  of  the  Solution 

In  the  previous  section  we  found  that  the  running  time  of  Algorithm  ECDFk  on  a  set  of  N 
points  in  k-space,  for  N  a  power  of  two,  is  given  by  the  expression 

w.N[.rwtn..i  «" 

In  this  section  we  shall  interpret  this  solution  to  see  what  it  tell  us  about  Algorithm  ECOFk. 
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In  many  applications  of  the  algorithm,  we  are  interested  in  the  case  that  k  Is  fixed  and  N 
grows;  we  can  then  use  the  approximation 

T(N,k)  -  N  ♦  KN  Ig^N).  (8) 

This  approximation  is  based  on  the  fact  that  (^)  -  p*Vk!  ♦  Wp^’*)  when  k  is  fixed  and  p 

grows.  It  is  interesting  to  remark  that  the  factor  1  /(k- 1 )!  in  the  expression  of  T(N,k)  was  still 
unknown,  and  that  it  explains  in  part  the  efficiency  of  MDC. 

The  above  analysis  holds  only  for  N  a  power  of  two.  The  running  time  of  Algorithm  ECDFk 
on  a  general  set  of  N  points  in  k-space  is  given  by  the  recurrence 

T(N,k)  -  T <r N/2l,k)  «•  T(lN/2J,k)  +T(N,k-l)  +  N  (9) 

T(N,2)  -  N  Ig  N 
T(l,k)  -  1. 

Note  that  when  N  is  a  power  of  two,  Equation  9  reduces  to  precisely  Equation  2  of  Section  3. 

We  shall  now  show  that  the  approximation  of  Equation  8  applies  also  to  this  general 

recurrence.  For  this,  we  define  the  first  and  second  differences  of  the  function  T  as 

aT(N,k)  -  T(N,k)  -  T(N-l.k-l) 

A2T(N,k)  -  AT(N,k)  -  AT(N-l.k-l). 

Using  the  recurrence  defining  the  function  T,  it  is  easy  to  find  the  recurrence  relations 
defining  the  differences.  It  becomes 

aT(N,k)  -  AT(N,k-l)  ♦  AT<rN/2l,k)  ♦  1 

A2T(N,k)  -  A2T(N,k-l)  for  even  N 

A2T(N,k)  -  A2T(N,k- 1 )  ♦  AT(fN/2l,k)  for  odd  N 

We  shall  prove  by  induction  on  k  and  N  that  aT(N,k)  is  positive  for  all  Nil  and  ki2.  For  k«2 
we  verify  this  easily  for  all  N  since  the  function  N  Ig  N  increases  with  N.  When  N-2,  we  have 

aT<2,k>  -  T«2,k>  -  T<  1  ,k)  -  T(l,k)  ♦  T(2,k-1)  ♦  1  >  0. 

Assume  that  the  first  difference  is  positive  in  dimension  k-1  for  any  N.  Then  a  sufficient 
condition  for  aT(N,k)  to  be  positive  is  that  aT(fN/2l,k)  is  positive,  and  so  on  until  we  need 

AT(2,k)  to  be  positive,  which  is  true.  We  have  therefore  proved  that  for  any  k  the  function 
T(N,k)  increases  with  N. 

Usino  fhe  same  argument  if  is  easy  fo  prove  that  the  second  difference  1$  always  positive, 
since  N  Ig  N  is  a  concave  function  of  N  and  A2T(2,k)  ■  a2T(2,2)  >0.  So  T(N^c)  Is  Itself  a 
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concave  function  of  N  for  any  fixed  k.  The  previous  results  will  be  sufficient  to  prove  that 
Equation  8  holds  for  all  N  and  for  any  fixed  k.  For  this  purpose,  we  shall  find  upper  and  lower 
bounds  of  T(N,k)  which  are  both  of  the  form  of  Equation  8.  Assuming  that  2P  <  N  <  2?*^,  and 
using  the  concavity  of  the  function  T(N,k),  we  may  bound  T(N,k)  by  S  and  I,  as  shown  in 

Figure  3.  Both  values  are  found  using  linear  interpolations  of  T  between  consecutive  powers 
of  two.  The  slopes  of  these  linear  interpolations  being  very  close,  we  find  for  S  and  ! 
expressions  which  are  equivalent  to  the  first  order  to 

D*'1 

s  -  I  -  N  ♦  #((N*2P)  p1"2) 

and  since  p  <  Ig  N  <  p*l,  the  expressions  of  S  and  I  are  identical  to  the  approximation  of 

T(N,k)  in  Equation  8.  Hence  we  have  proved  that  Equation  8  holds  for  any  N>1  and  kfc2. 


Figure  3s  Approximation  of  T(N,k)  using  the  concavity  of  the  function. 

This  allows  us  to  compare  Algorithm  ECDFk  with  the  naive  ECDF  algorithm  that  compares  all 
pairs  of  points,  whirh  we  raM  the  sequential  searching  ($$)  algorithm.  The  complexity  Of  this 
algorithm  is  0(k  N^).  For  fixed  k  and  large  enough  N  it  is  dear  that  ECDFk  is  better  than  SS, 
since  0(N  Ig^'^N)  <  0(N^).  In  some  applications,  however,  the  number  of  dimensions  Is  too 
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large  to  enable  us  to  use  tjie^apjjroximation  of  Equation  8.  We  then  approximate  the 
complexity  of  ECDFk  by  c  N  (  8  ^  )  and  the  complexity  of  $S  by  k  N^.  Figure  4  shows  in 
which  domains  each  algorithm  is  faster,  for  various  values  of  the  parameter  c.  Since  problems 
involving  more  than  280  (or  about  10^)  points  will  probably  never  be  processed,  Figure  4 
covers  the  actual  domain  of  values  for  k  and 


Figure  4:  Comparison  of  ECDFk  and  SS  algorithms. 

In  the  next  section  we  shall  analyze  other  MDC  algorithms  and  find  similar  expressions  for 
their  attributes.  The  same  interpretation  as  we  have  presented  for  this  problem  can  be  done 
for  the  problems  we  will  see,  and  will  yield  similar  conclusions. 

5.  Examples 

The  following  examples  analyze  the  preprocessing  time,  storage  and  query  time  of  two 
important  algorithms,  the  Maxima  Searching  and  the  ECDF  Searching  Problems.  A  description 
of  these  algorithms  may  be  found  in  Bentley  [1978]. 
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5.1.  Maxima  Searching 

The  first  problem  that  we  shall  analyze  is  that  of  the  maxima  searching  data  structure. 
There  are  three  attributes  of  this  data  structure  to  be  analyzed:  its  preprocessing  time  (how 
long  it  takes  to  build  the  structure),  its  query  time  (how  long  it  takes  to  search  in  the 
structure),  and  its  storage  (how  much  space  is  required  to  represent  the  structure).  The  first 
attribute  that  we  shall  analyze  is  the  preprocessing  time,  which  is  given  by  the  recurrence 

P(N,k)=2  P(N/2,k)+P(N,k-l)+0(N) 

P(N,2)=0(N) 

P(l,k)**0(l). 

We  can  transform  this  recurrence  using 
B(p,k)=P(2p,k)/2p  4  1 
into  the  following  system 

B(p,k)=B(p-l,k)4B(p,k-l) 

B(p,2)=2 

B(0,k)«2 

which  corresponds  to  a  net  similar  to  that  of  Figure  1  with  somewhat  different  boundary 
conditions.  The  solution  is 


(P*k-2\ 

B(p,k)  =  2  \  k-2  I 

and  we  find  the  preprocessing  time  to  be  equal  to 

fig  N  +(k-2)\ 

P(N,k)  -  2  N  l  k-2  I  -  N. 
For  fixed  k,  as  N  grows,  we  may  approximate  P(N,k)  by 

P(N,k)  «  2  N  4  #<N  lg*,-3N) 


We  now  turn  our  attention  to  the  storage  requirements,  which  are  defined  by 

S(N,k)-2  $(N/2,k)4S(N/2,k-l) 

S(N,2)«N 

S(l,k)-1. 

We  transform,  as  before,  by 
B(p,k)-S(2p,k) 
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B(p,k) 


X  2p*(M)  (k-i) 


2<i<k 

2  2-  (?) 


Osisk-2 


and  the  storage  is 


S(N,k)  -  N  X 
0<i<k-2 


We  may  remark  that  when  k  >  Ig  N  +2,  the  equation  becomes 

S(N,k)  =  X  2'  ('8i N)  -  3**  N  -  N1*  3 . 

0<i<k-2 

by  the  binomial  theorem.  An  intuitive  explanation  of  this  fact  is  that  for  k  large  enough,  we 
may  equate  S(N,k)  and  S(N,k-l)  and  write 

S(N,k)  =  2  S(N/2,k)  4  S(N/2,k-l)  *  3  S(n/2,k) 
whose  solution  is  precisely  S(N,k)  «  3*8  N  *  |\|ig  3 

At  the  opposite  end  of  the  spectrum,  for  fixed  k  and  increasing  N,  we  have 

S(N,k)  -  N  4  KN  lg‘‘-3N). 


The  last  attribute  of  the  maxima  searching  data  structure  Is  the  query  time.  Its  associated 
system  is 

CXN,k)«Q<N/2,k)4Q<N/2,k- 1  )40<  1 ) 

0(N,2)=lg  N 
QU,k)-l. 

The  tranformation 

B(p,k)-Q(2P,k)4l 

yields 

B(p,k)»B<p- 1  ,k)4B(p- 1  ,k- 1 ) 

B(p,2)-p4l 

B(0,k)-2. 

The  corresponding  net  is  similar  to  Figure  5,  but  there  is  only  one  way  for  each  diagonal 
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step.  So  we  find 

B(p,K)  -  (k-1 )  ♦  2  X  (?) 

OsiSk-2 

and  the  query  time  is 

fig  N  -1\  v  flgN\ 

Q(N,K)-l  K-1  1*2  2-  \  •  /-  1 

OSiSk-2 

For  fixed  K,  as  N  grows,  we  have 

CHN.K)  -  ♦  #(lgh'2N). 

We  must  remarK  that  the  query  time  admits  an  interpretation  similar  to  that  of  Section  4,  and 
similar  conclusion  as  for  the  comparison  with  the  obvious  SS  algorithm.  At  the  opposite,  the 
storage  and  the  preprocessing  time  are  0(N  K)  for  SS,  and  hence  better  than  those  of  the 
MDC  Maxima  Searching  algorithms. 

5.2.  ECDF  Searching 

The  second  problem  that  we  shall  study  is  the  ECDF  searching;  it  is  described  by  the  same 
three  attributes. 

The  preprocessing  time  and  the  storage  are  described  by  exactly  the  same  recurrences,  so 
we  shall  restrict  our  attention  to  the  preprocessing  given  by 

P(N,k)-2  P<N/2,k)*P<N/2,K-l)*WN) 

P(N,2)-N  Ig  N 
P(l,K)-l. 

The  function 

B(p,K)«P(2P,k)42P 
gives  us  the  new  system 

B(p,k)-2  B(p-l,k)4B(p-l,k-l) 

B(p,2)-(p4l)  2P 
B(0,k)-2 

whose  net  is  similar  to  Figure  5.  We  can  expand  the  net  by  Bfp.D^P4*  for  p>0.  The 
solution  is 
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B(p, 


,k)  -  2p'h*2  (k-1 )  ♦  2  S  2**'  (ft 

OsiSk-2 


The  preprocessing  time  and  storage  are  given  by 


Z  rM-i  j 

[  0£iSk-2  J 

Note  that  if  k  <  lg  N  +2,  then  P(N,k)  -  2  n'K  3  -  N.  For  fixed  k  and  large  N  the  following 
approximation  holds 


P(N,k)  -  N  ♦  UN  lgk-2N). 


The  last  parameter  to  be  analysed  is  the  query  time.  It  is  described  by  the  system 

Q(N,k)-Q(N/2,k)*Q<N/2,k- 1  > 

Q(N,2)-lg2  N 
Q<l,k)«l. 

We  use  the  new  function 
B(p,k)-Q(2P,k)*l 
and  we  have  to  solve 

B(p,k)=B(p-l,k)*Bfp-l,k-l) 

B(p,2)-pZ*» 

0(0, k  )-2. 

For  this  purpose,  we  expand  the  corresponding  net  by  noticing  that  B(p,l)  -  2p+l  and 
B(p,0)-2  are  convenient  for  p2l.  It  now  becomes 

0(p,k)  -  3  (k-l)  ♦  2  (Pk  )  ♦  2  X  (k-i)  . 

2<i£k 

-*(?).(?:!). *  Z  (?) 

OsiSk-2 


We  substitute  this  expression  in  the  initial  equation  to  find 

Wk>.2(vm\v).j  z  ru 

0£i£k-2 


Again,  we  may  remark  that  if  k  t  lg  N,  then  0(N,k)-2l*  N  -  N.  As  usual,  when  k  Is  fixed  and  N 


29  May  1979 


Combinatorial  Solutions  of  Recurrences 


-  14  - 


grows  we  find 


Q(N,k)  -  2  ^  ♦  #0gh‘,N). 


6.  More  General  Recurrences 

We  shall  now  study  how  we  can  solve  the  general  reduced  recurrence  of  the  form 
B(p,k)  -  a  B(p-l,k)  ♦  b  B(p-€,k-l) 

where  €-0  or  1.  The  values  of  B(p,k)  are  initially  definite  on  a  boundary  and  we  assume 
that  this  bound  u  of  the  form  of  that  in  Figure  6.  If  this  is  necessary,  we  just  consider  a 
subset  of  Eq  of  this  form. 


Figure  6:  The  sets  Eq  and  Ej  associated  to  a  net. 

We  now  define  the  set  Ef  to  be  the  set  of  points  which  are  not  in  Eq  and  which  can  be 
reached  from  Eq  in  one  step.  Then,  for  any  point  (p,k)  not  in  Eq,  the  function  B(p,k)  is  the 
sum  of  the  number  of  ways  from  Ej  to  (p,k),  each  way  weighted  by  the  value  of  B  at  its 


29  May  1979 


Combinatorial  Solutions  of  Recurrences 


-  15  - 


starting  point  rj-(pj,kj)  in  Ej.  So  we  find  the  general  solution 


B(p,k)  -  X  W  (p,k)  BfpA) 
r/E,  ' 


As  in  the  Section  3,  W(rs)(p,k)  denotes  the  number  of  ways  from  (r,s)  to  (p,k)  in  the  current 
net.  We  may  remark  that  the  sum  in  the  previous  equation  is  finite,  because  W(r>s)(p,k)  is  null 
out  of  a  finite  sub-net,  for  fixed  <p,k). 

It  may  happen  that  the  values  of  B  are  sufficiently  regular  on  Eq  for  allowing  us  to  expand 
the  net,  and  this  sometimes  reduces  to  a  sum  involving  only  a  constant  number  of  terms  in 
the  expression  of  B(p,k).^ 

Another  problem  is  now  to  reduce  the  recurrences  that  appear  in  MOC-problems  to  the 
previous  form.  The  general  paradigm  (divide  N  by  2  and  solve  the  problem  in  a  k-1  space) 
shows  that  the  natural  variables  are  Ig  N  and  k;  so  the  primary  change  of  variables  is  to  use 
p  m  |g  N.  We  therefore  solve 

B(p,k)  -  a  B(p-l,k)  +  b  B(p-«,k-l)  ♦  f(p,k) 

B(p,k)  **  known  on  a  bound. 

The  problem  is  to  find  a  particular  solution  A(p,k),  in  order  to  solve  the  reduced  equation 
verified  by  B-A.  The  only  way  seems  to  be  luck  and  trick,  and  no  general  method  can  be 
exhibited  here. 


7.  Conclusions 

In  this  section  we  will  briefly  review  the  contributions  of  this  paper.  One  of  the  main 
contributions  has  been  the  detailed  analysis  of  the  Multidimensional  Divide-and-Conquer 
algorithms  described  by  Bentley  [1978].  We  have  exhibited  precise  analyses  for  many  of 
those  algorithms,  accurate  to  within  an  (implementation-dependent)  constant  factor.  These 
analyses  show  that  the  algorithms  are  more  efficient  than  previously  thought  (it  was  not 
known  that  the  constant  of  proportionality  is  the  very  small  function  l/(k!)). 

In  addition  to  analyzing  particular  algorithms,  we  have  seen  a  set  of  general  tools 
applicable  to  the  analysis  of  algorithms  in  two  variables.  The  primary  analytical  tools  are  a 
set  of  useful  transforms  and  an  isomorphism  of  recurrences  and  path-counting  problems  on 
networks.  We  have  also  seen  a  number  of  tools  for  the  interpretation  of  recurrences;  these 
include  a  technique  for  showing  the  "smoothness”  of  the  resulting  function  (between  powers 
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of  two)  and  a  method  for  comparing  sophisticated  algorithms  with  more  straightforward 
solutions.  The  methods  that  we  have  seen  are  applicable  to  all  of  the  algorithms  described 
by  Bentley  [1978],  as  well  as  other  many  others  (such  as  in  Lee  and  Wong  [1979]). 
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